library(tigris)
library(acs)
library(leaflet)
library(mapview)
library(rgdal)
library(ggplot2)
library(htmltools)
library(maptools)
library(stringr)
api.key.install(key=mykey)createPercentMap <- function(title, table.name, variable.index) {
countylist <- c('17','21','25', '09')
shapefile <- tracts(state='MA', county=countylist, year=2018, class='sp')
# Span
myspan <- 5
# End year
myendyear <- 2018
mytable <- acs.lookup(endyear=2018, table.number=table.name)
# Limit the geography
# Convert the list of counties that is in text form to numeric form
countylist2 <- as.numeric(countylist)
# Define the geography to limit the data import
mygeo <- geo.make(state=25, county=countylist2, tract="*")
by.variable <- mytable[1]+mytable[variable.index]
mydata <- acs.fetch(endyear=myendyear, span=myspan, geography=mygeo, variable=by.variable)
mystate <- as.character(mydata@geography$state)
mycounty <- as.character(mydata@geography$county)
# add leading zeros to make the total number of characters 3
mycountypadded <- str_pad(mycounty, 3, pad="0")
mytract <- as.character(mydata@geography$tract)
acsgeoid <- paste0(mystate,mycountypadded,mytract)
mydatadf <- data.frame(acsgeoid, mydata@estimate)
colnames(mydatadf)=c("GEOID", "total", "variable")
mydatadf$variable.percent <- mydatadf$variable/mydatadf$total*100
mydatamerged <- geo_join(shapefile, mydatadf, "GEOID", "GEOID")
mydatamerged <- subset(mydatamerged, mydatamerged@data$variable.percent<100)
mypopup <- paste0("GEOID: ", mydatamerged$GEOID, "<br>", title, round(mydatamerged$variable.percent,1), "%")
mypal <- colorNumeric(
palette = "YlGnBu",
domain = NULL
)
myLAT <- 42.3013735
myLNG <- -71.1135916
myZOOM <- 14
myTILES <- "CartoDB.Positron"
mymap <- leaflet() %>%
addProviderTiles(myTILES) %>%
setView(myLNG, myLAT, zoom = myZOOM) %>%
addPolygons(data = mydatamerged,
fillColor = ~mypal(variable.percent), # this makes the fill color on the spectrum of the palette
color = "#b2aeae", # this is the color of the outline of the shapes
fillOpacity = 0.7, # how see through the shapes are
weight = 1,
smoothFactor = 0.2,
popup = mypopup) %>%
addLegend(pal = mypal,
values = mydatamerged$variable.percent,
position = "bottomright",
title = title,
labFormat = labelFormat(suffix = "%"),
bins=5)
return(mymap)
}
createFlatMap <- function(title, table.name, variable.index, prefix = "") {
countylist <- c('17','21','25', '09')
shapefile <- tracts(state='MA', county=countylist, year=2018, class='sp')
# Span
myspan <- 5
# End year
myendyear <- 2018
mytable <- acs.lookup(endyear=2018, table.number=table.name)
# Limit the geography
# Convert the list of counties that is in text form to numeric form
countylist2 <- as.numeric(countylist)
# Define the geography to limit the data import
mygeo <- geo.make(state=25, county=countylist2, tract="*")
by.variable <- mytable[variable.index]
mydata <- acs.fetch(endyear=myendyear, span=myspan, geography=mygeo, variable=by.variable)
mystate <- as.character(mydata@geography$state)
mycounty <- as.character(mydata@geography$county)
# add leading zeros to make the total number of characters 3
mycountypadded <- str_pad(mycounty, 3, pad="0")
mytract <- as.character(mydata@geography$tract)
acsgeoid <- paste0(mystate,mycountypadded,mytract)
mydatadf <- data.frame(acsgeoid, mydata@estimate)
colnames(mydatadf)=c("GEOID", "variable")
mydatamerged <- geo_join(shapefile, mydatadf, "GEOID", "GEOID")
mydatamerged <- subset(mydatamerged, mydatamerged@data$variable>0)
mypopup <- paste0("GEOID: ", mydatamerged$GEOID, "<br>", title, " ", prefix, mydatamerged$variable)
mypal <- colorNumeric(
palette = "YlGnBu",
domain = NULL
)
myLAT <- 42.3013735
myLNG <- -71.1135916
myZOOM <- 14
myTILES <- "CartoDB.Positron"
mymap <- leaflet() %>%
addProviderTiles(myTILES) %>%
setView(myLNG, myLAT, zoom = myZOOM) %>%
addPolygons(data = mydatamerged,
fillColor = ~mypal(variable), # this makes the fill color on the spectrum of the palette
color = "#b2aeae", # this is the color of the outline of the shapes
fillOpacity = 0.7, # how see through the shapes are
weight = 1,
smoothFactor = 0.2,
popup = mypopup) %>%
addLegend(pal = mypal,
values = mydatamerged$variable,
position = "bottomright",
title = title,
labFormat = labelFormat(prefix = prefix),
bins=5)
return(mymap)
}##
|
| | 0%
|
|= | 1%
|
|= | 2%
|
|== | 2%
|
|== | 3%
|
|== | 4%
|
|=== | 4%
|
|=== | 5%
|
|==== | 5%
|
|==== | 6%
|
|===== | 7%
|
|===== | 8%
|
|====== | 8%
|
|====== | 9%
|
|======= | 9%
|
|======= | 10%
|
|======== | 11%
|
|======== | 12%
|
|========= | 12%
|
|========= | 13%
|
|========== | 14%
|
|========== | 15%
|
|=========== | 15%
|
|=========== | 16%
|
|============ | 17%
|
|============ | 18%
|
|============= | 18%
|
|============= | 19%
|
|============== | 19%
|
|============== | 20%
|
|============== | 21%
|
|=============== | 21%
|
|=============== | 22%
|
|================ | 23%
|
|================ | 24%
|
|================= | 24%
|
|================= | 25%
|
|================== | 25%
|
|================== | 26%
|
|=================== | 27%
|
|=================== | 28%
|
|==================== | 28%
|
|==================== | 29%
|
|===================== | 29%
|
|===================== | 30%
|
|====================== | 31%
|
|====================== | 32%
|
|======================= | 32%
|
|======================= | 33%
|
|======================== | 34%
|
|======================== | 35%
|
|========================= | 35%
|
|========================= | 36%
|
|========================== | 37%
|
|========================== | 38%
|
|=========================== | 38%
|
|=========================== | 39%
|
|============================ | 39%
|
|============================ | 40%
|
|============================ | 41%
|
|============================= | 41%
|
|============================= | 42%
|
|============================== | 42%
|
|============================== | 43%
|
|=============================== | 44%
|
|=============================== | 45%
|
|================================ | 45%
|
|================================ | 46%
|
|================================= | 47%
|
|================================= | 48%
|
|================================== | 48%
|
|================================== | 49%
|
|=================================== | 49%
|
|=================================== | 50%
|
|=================================== | 51%
|
|==================================== | 51%
|
|==================================== | 52%
|
|===================================== | 52%
|
|===================================== | 53%
|
|===================================== | 54%
|
|====================================== | 54%
|
|====================================== | 55%
|
|======================================= | 55%
|
|======================================= | 56%
|
|======================================== | 57%
|
|======================================== | 58%
|
|========================================= | 58%
|
|========================================= | 59%
|
|========================================== | 59%
|
|========================================== | 60%
|
|=========================================== | 61%
|
|=========================================== | 62%
|
|============================================ | 62%
|
|============================================ | 63%
|
|============================================= | 64%
|
|============================================= | 65%
|
|============================================== | 65%
|
|============================================== | 66%
|
|=============================================== | 67%
|
|=============================================== | 68%
|
|================================================ | 68%
|
|================================================ | 69%
|
|================================================= | 69%
|
|================================================= | 70%
|
|================================================= | 71%
|
|================================================== | 71%
|
|================================================== | 72%
|
|=================================================== | 72%
|
|=================================================== | 73%
|
|=================================================== | 74%
|
|==================================================== | 74%
|
|==================================================== | 75%
|
|===================================================== | 75%
|
|===================================================== | 76%
|
|====================================================== | 77%
|
|======================================================= | 78%
|
|======================================================= | 79%
|
|======================================================== | 79%
|
|======================================================== | 80%
|
|======================================================== | 81%
|
|========================================================= | 81%
|
|========================================================= | 82%
|
|========================================================== | 82%
|
|========================================================== | 83%
|
|=========================================================== | 84%
|
|=========================================================== | 85%
|
|============================================================ | 85%
|
|============================================================ | 86%
|
|============================================================= | 87%
|
|============================================================= | 88%
|
|============================================================== | 88%
|
|============================================================== | 89%
|
|=============================================================== | 89%
|
|=============================================================== | 90%
|
|=============================================================== | 91%
|
|================================================================ | 91%
|
|================================================================ | 92%
|
|================================================================= | 93%
|
|================================================================== | 94%
|
|================================================================== | 95%
|
|=================================================================== | 95%
|
|=================================================================== | 96%
|
|==================================================================== | 96%
|
|==================================================================== | 97%
|
|===================================================================== | 98%
|
|===================================================================== | 99%
|
|======================================================================| 99%
|
|======================================================================| 100%# Create a map using leaflet and store that map in the variable mymap
subway <- readOGR(dsn="transitlayers/mbtarapidtransit", layer="MBTA_ARC") %>%
spTransform(CRS("+proj=longlat +ellps=GRS80"))## OGR data source with driver: ESRI Shapefile
## Source: "/Users/jea/Documents/summer_classes_2020/Rstuff/transitlayers/mbtarapidtransit", layer: "MBTA_ARC"
## with 138 features
## It has 4 fieldspal <- colorFactor(palette = c('blue','green','orange','red','grey'),
domain = c('BLUE','GREEN','ORANGE','RED','SILVER'))
# Define new MBTA icon
mbtaIcon <- makeIcon(
iconUrl = "https://upload.wikimedia.org/wikipedia/commons/thumb/6/64/MBTA.svg/240px-MBTA.svg.png",
iconWidth = 15, iconHeight = 15
)
# Make a purple circle icon for the commuter rail stations
crIcon <- makeIcon(
iconUrl = "https://upload.wikimedia.org/wikipedia/commons/thumb/8/80/Eo_circle_purple_circle.svg/240px-Eo_circle_purple_circle.svg.png", iconWidth = 15, iconHeight = 15)
# Load the MBTA stops shapefile
stops <- readOGR(dsn="transitlayers/mbtarapidtransit", layer="MBTA_NODE") %>%
spTransform(CRS("+proj=longlat +ellps=GRS80"))## OGR data source with driver: ESRI Shapefile
## Source: "/Users/jea/Documents/summer_classes_2020/Rstuff/transitlayers/mbtarapidtransit", layer: "MBTA_NODE"
## with 166 features
## It has 4 fieldscommRail <- readOGR(dsn="transitlayers/trains", layer="TRAINS_ARC") %>%
spTransform(CRS("+proj=longlat +ellps=GRS80")) %>%
subset(COMMRAIL == "Y", na.rm=TRUE)## OGR data source with driver: ESRI Shapefile
## Source: "/Users/jea/Documents/summer_classes_2020/Rstuff/transitlayers/trains", layer: "TRAINS_ARC"
## with 6595 features
## It has 28 fields
## Integer64 fields read as strings: TRACK LINE_CODE ASSET_ID# Import rail stations
# Subset to only include commuter rail stations
trainstops <- readOGR(dsn="transitlayers/trains", layer="TRAINS_NODE") %>%
spTransform(CRS("+proj=longlat +ellps=GRS80")) %>%
subset(C_RAILSTAT == "Y", na.rm=TRUE)## OGR data source with driver: ESRI Shapefile
## Source: "/Users/jea/Documents/summer_classes_2020/Rstuff/transitlayers/trains", layer: "TRAINS_NODE"
## with 387 features
## It has 6 fieldsbus <- readOGR(dsn="transitlayers/mbtabus", layer="MBTABUSROUTES_ARC") %>%
spTransform(CRS("+proj=longlat +ellps=GRS80"))## OGR data source with driver: ESRI Shapefile
## Source: "/Users/jea/Documents/summer_classes_2020/Rstuff/transitlayers/mbtabus", layer: "MBTABUSROUTES_ARC"
## with 923 features
## It has 11 fields
## Integer64 fields read as strings: CTPS_ROU_2# Subset the bus routes into "key routes" and "non-key-routes"
keybuslist <- as.character(c(1,15,22,23,28,32,39,57,66,71,73,77,111,116,117))
bus.key <- subset(bus, MBTA_ROUTE %in% keybuslist, na.rm=TRUE)
bus.nokey <- subset(bus, !(MBTA_ROUTE %in% keybuslist), na.rm=TRUE)
# Create the map
createPercentMap("Commute by Public transportation", "B08301", 10) %>%
addPolylines(data = bus.nokey, color="grey", weight=1) %>%
addPolylines(data = bus.key, color="yellow", weight=2) %>%
addPolylines(data = commRail, color="purple", weight=2) %>%
addPolylines(data = subway, color=~pal(LINE)) %>%
addMarkers(data = stops, icon=mbtaIcon, label = ~htmlEscape(STATION), labelOptions = labelOptions(noHide = T)) %>%
addMarkers(data = trainstops, icon=crIcon, label = ~htmlEscape(STATION), labelOptions = labelOptions(noHide = T))